#ifndef MATRIX_H  
#define MATRIX_H  
  
#include <stdio.h>  
#include <iostream>  
#include <vector>  
using namespace std;  
  
/**  
 * @brief 模板类Matrix用于表示矩阵并提供相关操作。  
 *   
 * @tparam T 矩阵元素的数据类型  
 */  
template<typename T>  
class Matrix {  
	private:  
		int n, m;  
		vector<vector<T>> data;  
  
	public:  
		/**  
		 * @brief 输入矩阵的元素。  
		 */  
		void input();  
  
		/**  
		 * @brief 输出矩阵。  
		 */  
		void print();  
  
		/**  
		 * @brief 将矩阵设置为单位矩阵。  
		 *   
		 * @note 仅适用于方阵。  
		 */  
		void setIdentity();  
  
		/**  
		 * @brief 默认构造函数。  
		 */  
		Matrix() {};  
  
		/**  
		 * @brief 构造一个指定大小的矩阵，矩阵元素赋初值为0。  
		 *   
		 * @param n 矩阵的行数  
		 * @param m 矩阵的列数  
		 */  
		Matrix(int n, int m);  
  
		/**  
		 * @brief 析构函数。
		 */  
		~Matrix() {};  
  
		/**  
		 * @brief 重载"+"运算符，实现矩阵的加法。  
		 *   
		 * @param other 另一个矩阵  
		 * @return Matrix<T> 两个矩阵相加的结果  
		 */  
		Matrix operator+(const Matrix& other) const;  
  
		/**  
		 * @brief 重载"*"运算符，实现矩阵的乘法。  
		 *   
		 * @param other 另一个矩阵  
		 * @return Matrix<T> 两个矩阵相乘的结果  
		 */  
		Matrix operator*(const Matrix& other) const;  
  
		/**  
		 * @brief 重载"*"运算符，实现矩阵与标量的乘法。  
		 *   
		 * @param num 标量值  
		 * @return Matrix<T> 矩阵与标量相乘的结果  
		 */  
		Matrix operator*(const T num) const;  
  
		/**  
		 * @brief 重载"*"运算符，实现标量与矩阵的乘法。  
		 *   
		 * @param num 标量值  
		 * @param A 矩阵  
		 * @return Matrix<T> 标量与矩阵相乘的结果  
		 */  
		template<typename S>  
		friend Matrix<S> operator*(const double num, const Matrix<S>& A);  
  
		/**  
		 * @brief 计算矩阵的乘方。  
		 *   
		 * @param p 指数  
		 * @return Matrix<T> 矩阵的乘方结果  
		 */  
		Matrix Pow(int p) const;  
  
		/**  
		 * @brief 计算矩阵的逆矩阵（未具体实现）。  
		 *   
		 * @note 仅适用于方阵。  
		 * @return Matrix<T> 矩阵的逆矩阵  
		 */  
		Matrix inverse();  
};  
  
template<typename T>  
void Matrix<T>::input() {  
	for (int i = 0; i < n; ++i) {  
		for (int j = 0; j < m; ++j) {  
			cin>>data[i][j];  
		}  
	}  
}  
  
template<typename T>  
void Matrix<T>::print() {  
	for (int i = 0; i < n; ++i) {  
		for (int j = 0; j < m; ++j) {  
			cout<<data[i][j]<<" ";  
		}  
		cout<<endl; 
	}  
}  
  
template<typename T>  
void Matrix<T>::setIdentity() {  
	if (n != m) printf("Error: matrix sizes mismatch!");  
	for (int i = 0; i < n; ++i) {  
		for (int j = 0; j < m; ++j) {  
			data[i][j] = (i == j);  
		}  
	}  
}  
  
template<typename T>  
Matrix<T>::Matrix(int n, int m) : n(n), m(m) {  
	data.resize(n, vector<T>(m, 0));  
}  
  
template<typename T>  
Matrix<T> Matrix<T>::operator+(const Matrix& other) const {  
	if (n != other.n || m != other.m) {  
		printf("Error: matrix sizes mismatch!");  
		return Matrix<T>(0, 0);  
	}  
	Matrix<T> result(n, m);  
	for (int i = 0; i < n; ++i) {  
		for (int j = 0; j < m; ++j) {  
			result.data[i][j] = data[i][j] + other.data[i][j];  
		}  
	}  
	return result;  
}  
  
template<typename T>  
Matrix<T> Matrix<T>::operator*(const Matrix& other) const {  
	if (m != other.n) {  
		printf("Error: matrix sizes mismatch!");  
		return Matrix<T>(0, 0);  
	}  
	int r = other.m;  
	Matrix<T> result(n, r);  
	for (int i = 0; i < n; ++i) {  
		for (int j = 0; j < r; ++j) {  
			for (int k = 0; k < m; ++k) {  
				result.data[i][j] += data[i][k] * other.data[k][j];  
			}  
		}  
	}  
	return result;  
}  
  
template<typename T>  
Matrix<T> Matrix<T>::operator*(const T num) const {  
	Matrix<T> result(n, m);  
	for (int i = 0; i < n; ++i) {  
		for (int j = 0; j < m; ++j) {  
			result.data[i][j] = data[i][j] * num;  
		}  
	}  
	return result;  
}  
  
template<typename T>  
Matrix<T> operator*(const double num, const Matrix<T>& A) {  
	Matrix<T> result(A.n, A.m);  
	for (int i = 0; i < A.n; ++i) {  
		for (int j = 0; j < A.m; ++j) {  
			result.data[i][j] = num * A.data[i][j];  
		}  
	}  
	return result;  
}  
  
template<typename T>  
Matrix<T> Matrix<T>::Pow(int p) const {  
	if (n != m) {  
		printf("Error: This matrix is not square!");  
		return Matrix<T>(0, 0);  
	}  
	if (p == 0) {  
		Matrix<T> result(n, n);  
		result.setIdentity();  
		return result;  
	}  
	if (p == 1) return *this;  
	Matrix<T> result = Pow(p >> 1);  
	result = result * result;  
	if (p & 1) result = result * (*this);  
	return result;  
}  
  
#endif  